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Abstract 

The observation of astrophysical neutrinos requires a detailed understanding of the atmo- 
spheric neutrino background. Since neutrinos are produced in meson decays together with a 
charged lepton, important constraints on this background can be obtained from the measure- 
ment of the atmospheric muon flux. Muons, however, can also be produced as u + /i~ pairs 
by purely electromagnetic processes. We use the Z-moment method to study and compare 
the contributions to the atmospheric muon and neutrino fluxes from different sources (n/K 
decay, charmed and unflavored hadron decay, and photon conversion into a muon pair). We 
pay special attention to the contribution from unflavored mesons (w, n', p° , w and <f)). These 
mesons are abundant in air showers, their lifetimes are much shorter than those of charged 
pions or kaons, and they have decay branching ratios of order 1CP 4 into final states contain- 
ing a muon pair. We show that they may be the dominant source of muons at > 10 3 
TeV. 



1 Introduction 



Atmospheric muons and neutrinos are produced in the showers of high energy cosmic rays in the 
Earth's atmosphere. Neutrinos can change flavor during their travel from the creation to the 
detection points, and therefore the observation of their fluxes allows to study their masses, mix- 
ings and interactions. Above a minimum energy of a few GeV muons reach the ground, where 
they can also be observed. Measurements of the atmospheric muon flux provide information 
about the primary cosmic ray fluxes and about the properties of high-energy hadronic inter- 
actions. Moreover, the measurements of the atmospheric muons can be used to obtain precise 
estimates of the fluxes of and before flavor oscillation effects. In fact, these measurements 
have played an important role for the interpretation of the data and the determination of the 
flavor oscillation parameters. The dominant source of muons and muon neutrinos is the weak 
decay of charged pions and kaons. These decays always produce £ + f£ or £~U£ pairs, implying a 
correlation between charged-lepton and neutrino fluxes that can be robustly predicted. 

An important goal of present and future experiments is the detection of high-energy neutrinos 
produced in galactic and extragalactic astrophysical sources [U [2] . Neutrinos propagate without 
significant losses (excluding the cosmological redshift) from very distant sources, and one of 
largest expected signals is an isotropic diffuse flux generated by the ensemble of all extragalactic 
sources in the universe. For neutrino astronomy, therefore, the atmospheric neutrino fluxes 
constitute a disturbing background that must be subtracted from the signal. Measurements of 
the muon flux can help in the determination of this background. 

At very high energy it is possible, and indeed virtually certain, that the decay of charged 
pions and kaons does not remain the dominant source of atmospheric muons and neutrinos. 
The reason is that these particles are relatively long lived and, because of the Lorentz time 
dilatation, their decay probability at high energy is strongly suppressed. As an illustration, the 
decay length of a 10 TeV charged pion is more that 500 Km, around 100 times larger than its 
interaction length in air. The contribution of particles with a shorter lifetime is therefore likely 
to become dominant. In particular, the contribution of charmed hadrons is a natural candidate. 
These particles have large (~ 10%) branching ratios into semileptonic modes and a lifetime 
r ~ 10 -12 s, implying a decay probability of order 1 up to energies around 10 7 GeV. 

Recently, however, it has been suggested [3] that the dominant source of atmospheric muons 
of very large energy could be the electromagnetic decay of unflavored mesons into [i^ fj,~ pairs. 
These mesons are particles of type q/Qf, constituted by a quark and an anti-quark of the same 
flavor. Neglecting heavy quarks one has 3 scalar (tt°, t\ and rf) and three vector uj and eft) 
mesons of this kind. They decay mostly into pions and photons via strong or electromagnetic 
interactions with a very short lifetime. All of them except the neutral pion (which is below 
threshold) can decay into a pair, sometimes together with a photon or a neutral pion 

(77 —7- ^ + /x _ 7 or uj —7- fi + fi~ir°), with small branching ratios of order 10 -4 . These rare decay 
modes have been neglected in essentially all calculations of atmospheric muons. However, if 
the multiplicities and the energy spectra of charged pions and unflavored mesons are roughly 
similar, then they will become the dominant source of atmospheric muons when the average 
decay probability of charged pions is suppressed by a factor (~ 10~ 4 ) of the same order as the 
branching ratios into the rare modes that contain muons. 

In this article we do a critical review of the different contributions to the lepton fluxes at high 
energies. In Section 3 we estimate the conventional lepton fluxes using the so called Z-moment 
method, that provides simple analytic expressions. Then we focus on the contribution to the 
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muon flux from the decay of unflavored mesons. In sections 5 and 6 we evaluate, respectively, 
the contribution from the prompt decay of charmed hadrons and from 7 conversion into a 
pair. Finally, we summarize the uncertainties and the implications of our results. 

2 Components of the atmospheric lepton fluxes 

The atmospheric flux of a lepton type £ can be described as the sum of contributions from the 
decay of different unstable particles: 



In this equation E is the energy and 6 the zenith angle of the lepton, and the summation 
runs over all possible parent particles. In this work with lepton we refer only to neutrinos and 
muons. The production of electrons and positrons is dominated by photon conversion in the 
electromagnetic field of the air nuclei (7Z — > e + e~ Z, with Z the electric charge of the nucleus) 
and will not be discussed here. 

The parent particles that are the source of atmospheric leptons can be naturally divided 
in three classes. The first class (standard contribution) includes charged pions and kaons that 
decay via charged-current weak interactions into lepton pairs: (e + ^ e ), (/u + f M ) and the charge 
conjugate states. The observed atmospheric fluxes can at the present time be entirely attributed 
to this standard contribution. 

The second source of atmospheric leptons (charm contribution) is the weak decay of particles 
that contain a charm (anti)-quark. These decays also generate leptons in {ev e ) and (/jlu^) 
pairs. In addition, Dg mesons (scalar mesons with a cs or sc quark content) have a decay 
branching ratio of ~ 6.4% into the 2-body mode t + u t (t~V t ); the subsequent decay of the 
t lepton generates a second tau (anti)-neutrino. This chain decay process is the main source 
of atmospheric u T and V T . The contribution of charmed particles to the atmospheric lepton 
fluxes is subdominant and currently undetected. It is however expected that this mechanism 
will overtake the standard contribution at sufficiently high energy. 

A third class of parent particles can contribute to the flux of atmospheric muons (but not of 
neutrinos). This unflavored contribution is due to the decay of the unflavored mesons 77, n' , p° , 
lu and 4>. These particles have small (order 10 -4 ) branching ratios into final states that include 
a pair. The possible significance of this contribution has been discussed in [3], and it will 

be critically analyzed in the following. 

Muon pairs can also be directly produced in Drell-Yan processes and in photon conversions 
of type 7Z — > fi + [i~Z. This last process, despite being suppressed by a factor (m e /m^) 2 ~ 
2.3 x 10 -5 with respect to the production of e + e~ pairs, is potentially interesting. Its contribution 
to the muon atmospheric flux will be indicated as cpy! and discussed later. 

Taking into account these different sources, the muon and neutrino fluxes can then be ex- 
pressed as the sum of four and two components, respectively: 




(1) 
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(2) 



M E > ) 
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3 The standard contribution 



Simple analytic expressions for the atmospheric lepton fluxes produced by the decay of charged 
pions and kaons are described in the textbook [3] by Gaisser (see [5] for additional details). 
These expressions are obtained under three simplifying assumptions: 

(i) The interaction lengths of all hadrons (labeled by k) are taken constant, neglecting 
their energy dependence. 

(ii) The inclusive spectra of secondary particles j created by the projectile particle in a 
hadronic interaction with an air nucleus satisfy the scaling condition: 

where x = Ej/E^, is the total inelastic cross section and Fj t j{x) is the number density 
of particles j carrying a fraction x of the initial energy after the collision. 

(iii) The primary nucleon fluxes are simple power laws of exponent a: 

cf> p (E ) =K po E^ a ; 

(5) 

MEo) =K n E^ a = K (1 - Po ) E^ a . 

In the low-energy limit, when the parent particle (a charged pion or kaon) decays with 

probability close to one, the lepton fluxes are isotropic and have a power-law energy spectrum 

(i) 

with the same exponent a as the primary nucleon fluxes. The component 4>f takes then the 
form 

4>f{E, 6) = K E~ a Aj(a) Z j£ {a) . (6) 

Aj(a) is the ratio between the number of nucleons that reach the Earth with energy in the 
interval (E,E + dE) and the number of particles of type j produced in the same energy interval 
by primary or secondary particles. This quantity is less than one even if primary nucleons tend 
to generate many secondary particles of type j, because these particles are produced with lower 
energy while the ratio is performed at a fixed E. The second quantity in Eq. ([6]), Zj£(a), is 
analogous, it relates the lepton flux with the flux of its parent particle j. It can be calculated 
as the (a — l)-moment of the inclusive spectrum Fji(x) of lepton £ from the decay of j: 

Z ji (a)= f'dxx - 1 F je (x), (7) 
J o 

where x = E^/Ej and Fji(x) is taken in any frame where the parent particle is ultrarelativistic. 
The quantity Aj{a) includes proton and neutron contributions: 

Aj (a) = po A pj (a) + n A nj (a) . (8) 

It is straightforward (see [HE]) to obtain these contributions in terms of Z-factors: 

A ■ + A ■ - Zpj ^ Znj (Q) 
r ^pp ~r ^pn 
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where the dependence on a is implicit. Again, the Z-factor Zkj(ct) is just the (a — l)-moment 
of F kj (x), 

Z kj (a) = [ dx x a - 1 F kj {x). (10) 



At high energy the decay probability of pions and kaons is suppressed because of the Lorentz 
time dilatation. When the decay probability of the parent particle is small, using the assumptions 
(i), (ii) and (hi) it is possible to express the lepton flux from j-decay as 

^\E,6) = K E' a £ii?zc g thW Bj(a) Z jt (a + 1) . (11) 

This energy spectrum is also a power law, but its slope is a unit steeper than in the primary 
nucleon flux. In addition, it has the strong dependence on the zenith angle described by -F zen ith(0) 
(shown in Fig. [TJ. For 9 < 60° this function is well approximated by a "secant law", 

^zenithW - — a , (12) 
COS V 

whereas for larger zenith angles it keeps growing monotonically, reaching at 6 ~ 90° a value 
close to 10. To a good approximation the zenith angle dependence obtained in this high energy 
limit is universal, it does not depend on the parent particle type or the details of the hadronic 
interactions. The quantity Ej is the critical energy for particle j: 

,,J^ (13) 

that corresponds to the condition where the decay length is equal to the scale height of the air 
density in the stratosphere h^ ~ 6.36 Km (p(h) cx e~ h ^ h °). The critical energies for ir^, Kl and 
are approximately 115, 210 and 850 GeV, respectively. After averaging over the creation 
position, the decay probability of a particle of type j energy E and zenith angle 6, for large 
energy takes the asymptotic form: 

p _ gj -^zenith (0) Bj(a) £j F zenith (9) 

dec " FT ' Aj(a) ~ ~ 1 j 



Note also that the decay Z-factor in Eq. (11) is calculated for the argument (a + 1). Finally, 
Bj(a) is analogous to the quantity Aj(a) defined in the low-energy limit but includes effects 
due to the distribution of the creation point of particle j. It can be separated as 

Bj(a) =p B pj (a) +n B nj (a) . (15) 

Including only the particles produced in nucleon interactions one obtains: 

where the dependence on a is implicit and 

A " = i z\z ■ (17) 

1 ^pp ~r ^pn 
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For pions, the inclusion of regeneration effects (the contribution of pions produced in pion 
interactions) yields the result: 



B pn ± (a) 



±B r . 



: (a) 



Pj Znj 



1 — Zpp =p z. 



J pn 



A± 



log 



At 



(18) 



with: 



A 



± 



1-Zt 



(19) 



Analogous expressions for kaons that include the effect of regeneration are discussed in [5]. 



Eq. (11) has been obtained under the hypothesis that the decay probability of the parent 

j Zenith (&)■ Correspondingly, the range of validity of the flux 
A useful expression that interpolates between the asymptotic 



particle j is small, i.e., E 3> e 
in Eq. ^ is E < Ej F zenitl 
fluxes in ^ and (11) is: 



(K E~ a ) Aj(a) Z jt (a) 



1 + 



E 



Aj(a) Zji(a) 



^zenith (0) Bj{a) Z jt (a+1) 



(20) 



The analytic expressions described above have a limited validity, since assumptions (i), (ii) 
and (iii) are not rigorously correct. The hadronic cross sections grow slowly with energy; the 
energy scaling Q is not exact; and the high-energy cosmic ray flux cannot be described as a 
simple power law because of the steepening at the cosmic ray knee (-Ekncc — 3 x 10 6 GeV). In 
a first approximation, the energy dependence of hadronic interaction lengths and the violation 
of scaling in inclusive particle distributions can be taken into account just by considering the 
quantities and Z^(a) as (slowly varying) functions of the lepton energy. The calculation of 



the lepton fluxes for an arbitrary shape of the primary flux is discussed in Subsection 3.2 



3.1 Numerical estimate 



The objective in this work is to compare the different contributions to the lepton fluxes at very 
high energy (Eg > 10 TeV). To estimate the standard contribution from pion and kaon decays 



one can therefore use Eq. (11). Summing over all parent particles this contribution (for zenith 
angles 9 < 60°) takes the form: 



/.stand 



(E,t 



E £ (a) 



(K E~ a ) E cos 6' 
where the constant Eg(a) has dimension of energy and is given by 

E e (a)= e 3 B M Zji{a+1). 



(21) 



(22) 



For a numerical estimate of the quantities Eg(a), we first have obtained the inclusive particle 
spectra from a MonteCarlo simulation generated by the code Sibyll [6], and then we have deduced 
the hadronic Z-factors from numerical integrations. The interaction lengths in air have been 
taken from the PDG fit [7] of hadron-nucleon cross sections, and we have used a Glauber 
formalism [8] to compute the cross sectios on a nuclear target. The hadronic Z-factors at 
Eg ~ 10 6 GeV that we obtain are shown in Table [l] The characteristic energy for muons, 
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E, [GeV] 


a = 2.7 


a = 3.0 


(h + + ij<~) 


9.3 


4.9 




3.4 


1.7 


(u e + Z7 e ) 


0.17 


0.10 


M + 


5.2 


2.8 


/i 


4.1 


2.1 




2.2 


1.1 




1.2 


0.6 




0.10 


0.06 




0.07 


0.04 



Table 1: Characteristic energies E^ (see definition in Eq. (21)) for the lepton fluxes calculated 



using the Sibyll interaction model. The quantities are evaluated for Eg ~ 10 GeV. 
summing over both charges, is 

E„(a = 2.7) ~ 9.3 GeV , 

(23) 

E„(a = 3.0) ~ 4.9 GeV . 

These estimates are given for the value of a relevant below and above the cosmic-ray knee. 

The energy dependence of the inclusive particle spectra in the fragmentation region (that 
dominates the integrand of the Z-factor expression) in the Sibyll code is weak, and the char- 
acteristic energies for Ei ~ 10 6 GeV (given in Table [l]) are approximately 25% lower than the 
results obtained at Ei ~ 10 GeV, where the calculation is in reasonable agreement with ob- 
servations. Because of our poor understanding of hadronic interactions one should attribute 
around ~ 30% of systematic uncertainty to the numerical estimates in Table [I] and Eq. ( 23 ) , 
but a larger error cannot be excluded. 

Finally we note that the absolute prediction for the lepton fluxes also depends on the estimate 
of the primary nucleon flux j9j HUJ [H] . This introduces perhaps the largest uncertainty, because 
the cosmic-ray composition above the knee is poorly known. 

3.2 Lepton yields 

To compute the lepton fluxes when the primary spectrum is not a power law it is necessary to 
introduce the lepton yields. The lepton yield Y p r n \_y£(E, Eq, 9) of a primary nucleon of energy Eq 
and zenith angle 9 gives the average number of leptons £ of energy E observable at ground level 
(where to a good approximation the shower has completed its development) per unit energy: 

Y p[n) ^(E, Eo, 9) = ( £, £ 0> e) . (24) 

The lepton flux can then be calculated with a simple integration: 



CO 



(E, 9) = / dEo [<P p (Eo) Yj^ t (E, E Q , 9) + <f> n (E ) Y n ^ e (E, E , 9)] . (25) 



E 



Y p ( n ^£ can be written as a sum of terms associated to the production and decay of different 
parent particles: 

Y p{n) ^{E, E , 9) = J2 Y$i)-JE, E , 9) , (26) 
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where the index j runs over all species with decay modes containing the lepton t. 

Using assumptions (i) and (ii) these yields take simple scaling forms. In the low-energy limit 
described in the previous section one has 



Eq, i 



1 



Q%(x) = ^- [G Nj ®F jt ]{x) 



1 



E 
1 

Eq Jo 



Eq 



dx\ I dx 2 Gnj(xi) Fje(x 2 ) S(x - x\ x 2 ) , 



(27) 



with x = E/Eq. Gpjj(x) dx is the average number of particles j created in the shower (in 
interactions of N and secondary particles with air nuclei) with a fraction of energy in the 
interval [x,x + dx]. The function G^j(x) is related to A^j{a) by the expression: 



A Nj (a) 



dx x a 1 Gnj{x) 



(28) 



GNj{x) can therefore be calculated as the inverse Mellin transform of A^j{a) 
In the high-energy limit the lepton yield takes the form 



yJ\i{E, E , 



E% cos 6 



Eq COS I 



Ni 

1 



El cos 9 



H 



dx\ 







dx-. 



X\ 

H Nj (xi) 

Xl 



(x) 



Fji{x 2 ) 5(x-x 1 x 2 ) 



(29) 



Again, H^j{x) includes the production of hadrons j in primary and secondary interactions inside 
the shower started by N. This function is related to B^j{a) by 



B 



Nj 



dx x 



a-l 



H 







Nj{ 



(30) 



and can therefore be calculated as its inverse Mellin transform. 

It is straightforward to ckeck that for a power-law primary flux, using Eqs. (27) and (29) to 
perform the integration in Eq. (25), one recovers for the lepton fluxes the results in Eqs. ([6]) and 

dub. 



4 Atmospheric muons from unflavored mesons 

The decay of the unflavored mesons {rj, r/, p°, u, (j)} into final states that contain a fi + 
pair also contributes to the atmospheric muon fluxes. It is straightforward to use the methods 
outlined in the previous section to estimate this contribution. 

Using the assumptions (i), (ii) and (iii) the muon flux from rj mesons is 



(E) 



(K E~ a ) 



A n (a) Z m (a) 



+ 



1 — Znn (1 — Znn) (1 — Zj, 



(31) 



with similar expressions for all other unflavored mesons. Eq. (31 ) provides an estimate of the sum 
of the /i + and [i~ fluxes, taking into account the production of eta mesons both in nucleon and 
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pion interactions (we neglect the smaller contribution from kaon interactions). The contributions 
from the decay of unnavored mesons to the positive and negative muon fluxes are identical. 

The decay Z-factors are computed from the measured branching fractions into states with 
muon pairs and from the shape of the muon energy spectra (in Fig. [2]). Numerical estimates for 
these factors are given in Table [2j These estimates have an uncertainty of order 20% associated 
to the experimental error in the measurement of the relevant branching ratios. 

The calculation of the Zn v and Z nr) factors obviously requires the modeling of unnavored 
meson production in nucleon and (less critically) pion interactions. The hadronic Z-factors 



entering Eq. (31) are 

Znt] = Zprj = Z n >q , 

Z-nri = Ztt+jj = Z n - n , 

Znn = Zpp + Zpn = Z nn + Z np , 

where we have left the dependence on a implicit and have used isospin symmetry. 

The inclusive meson spectra have been obtained from a Sibyll |6J Montecarlo simulation (see 
Fig. [3]) , and we have then evaluated the corresponding Z-factors through numerical integration. 
The results for a = 2.7 and a = 3 are listed in Table |3j The combination 

Zunflav(a) = Z Nj (a) (33) 

is also plotted as a function of a in Fig. |4j 

The contribution to the muon flux from unnavored mesons results from the addition 

^3^= J2 Ma)Z 3i (a) = C^(a). (34) 

je{r?,r/,p°,w,0} 

From the Z-factors in Tables [2] and |3] we obtain 

c unfiav( Q = 2 _ g 2 x 1Q -6 ^ 

(35) 

C^ nfiav (a = 3.0) ~ 3.1 x 10~ 6 . 

In the regions of the spectrum where the primary nucleon flux is not well described by a 
single power law it is possible to compute 0^ nfl (.E) from the muon yields. The results of such 
calculation that we obtain using the Sibyll Montecarlo code are given in Fig. [5] The top line 
there shows the all-nucleon primary flux. The thick red line is our estimate of the muon flux 



X 


Z^(2.7) (xlO- 4 ) 


Z^(3.0) (xl0~ 4 ) 


V 


1.37 


1.12 


rf 


0.43 


0.35 


P° 


0.33 


0.30 




1.00 


0.86 





2.15 


1.93 



Table 2: Decay Z factors for unnavored mesons. 
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Z p j(2) 


Z nj {2) 


Z pj (2.7) 


Z vj {2.7) 


Z pj (3.Q) 


Z WJ {3.0) 


V 


0.066 


0.094 


0.014 


0.029 


0.0087 


0.021 




0.052 


0.074 


0.013 


0.027 


0.0086 


0.020 


p° 


0.054 


0.077 


0.013 


0.026 


0.0082 


0.019 




0.040 


0.060 


0.010 


0.021 


0.0066 


0.016 




0.0019 


0.0020 


0.00038 


0.00047 


0.00022 


0.00029 


All 


0.22 


0.31 


0.051 


0.10 


0.032 


0.076 



Table 3: Z factors for the production of unfiavored mesons in proton and pion collisions with 
an air nucleus. The inclusive spectra are calculated with the Sibyll Montecarlo code [6j 



from the electromagnetic decay of unfiavored mesons. We include the conventional /i + + 

v a + T^fi and v e + v e fluxes (from the vertical direction) taking into account only the decay of 

charged pions and kaons. The data points are from the L3 detector |12j . 

Inspection of Fig.[5]shows that for the vertical direction {0 = 0) the contribution of unfiavored 
muons overtakes the standard contribution from pion and kaon decay at ~ 1.6 x 10 6 GeV. 
This result can be also estimated combining Eqs. (23) and (35). The unfiavored contribution is 
isotropic, while the standard contribution grows with increasing zenith angle proportionally to 
^zenith (#)• 



4.1 Unfiavored meson production 

We would like to analyze the uncertainties in the modeling of unflavored-meson production in 
hadronic interactions, and in this subsection we give a qualitative discussion. 

The most important quantities for the prediction of the muon flux are the fraction of the 
projectile energy carried by the unfiavored mesons produced in the collision and the shape of 
their energy spectra. To estimate the energy fraction (xj) = (Ej)/Eo carried by the meson type 
j, one can make some simple considerations. The initial energy Eq of the projectile is divided 
among three classes of particles: baryons, antibaryons and mesons, 

E ^ E qqq + E Wj + E q g . (36) 

The energy E q g carried by mesons is in turn subdivided among different particle types: 

Egg= Y, ^ 

j'£{mesons} 

The summation is over the "primary mesons (before the decay of unstable particles) . It is a good 
approximation to neglect heavy quarks and include in the sum only the 18 mesons that compose 
the scalar and vector nonets of SU(3). Meson production can be modeled as a two-step process: 
in the first step pairs are created and recombined between each other and with the valence 
quarks of the interacting nucleons; in the second step the states qj~q k are projected into physical 
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mesons, whose quark content is known. The scalar unflavored mesons have the quark content 

7r° = —= (uu — dd\ ; 
v2 

rj = - (uu + dd) - -j= ss ; (38) 

rf = - (uu + dd) H — -= ss , 
2 \/2 



whereas for the vector mesons: 



1 
1 



uu — dd) ; 

dd) ; (39) 



ss . 



Assuming that each light flavor combination qjq k has on average the same energy and ne- 
glecting meson mass differences, the fraction of energy carried by the different species can be 
calculated in terms of just two parameters: the probability P s of producing an ss pair (with 
Pu = Pd = (1 — -fs)/2), and the probability -P S calar to project the state qjq k into a spin-0 meson. 
The energy fractions can then be estimated as 



(E n ) (E v 



E qq E qq 



P. 



scalar 



8 2 



-,(!-*_,.) <i^£ (40) 



E - R - 



(1 -^scalar) -P< 



2 



For completeness, the average energy taken by the other mesons is 



(E n o) ^ (E n+ ) ^ (E n -) ^ (1-P t 



\2 



-^scalar 



Eqq E qq E qq i 

(E k+ )„(Ek-) „,{Ek l ) „,{E Ks ) ^ (1-P S )P S 

TP _ TP — TP — TP — scalar n 

^qq QQ QQ QQ 



For the corresponding vector particles (p and K*) one can use Eq. (41) with the substitution 

-^scalar ^ 1 P, scalar- 

For a numerical estimate we can use E qq / Eq ~ 0.6 (with most of the remaining energy carried 
by one leading baryon) and, following Field and Feynman (13] and the Lund fragmentation 
algorithm |14j . P sca i ar — 0.5 and P s ~ 0.13. With these assumptions the 5 unflavored mesons 
carry together an energy fraction (-E U nflav) /Eq ~ 0.18. This estimate depends only weakly on the 
values chosen for P S caiar and P s . The minimum value, (-E U nflav) — 0.13, is obtained for P S caiar = 1 
and P s = 0. 

It is reasonable to expect that the energy spectrum of unflavored mesons is similar or slightly 
harder than the spectrum observed for pions (which may come from a longer decay chain of 
unstable primary mesons) . A simple 2-parameter form for the inclusive energy spectrum F^j (x) 
is: 

F Nj (x) = ( Xj ) (1 + nj) { X ->- . (42) 
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where (xj) is the energy fraction carried by particle j and rij a shape parameter in the range 



3-4. The Zjvj moments corresponding to (42) are 



*»,(.)-<„) £±«jW , (43 ) 

T(n i + a) 

Note that Zj-j(2) = (xj) and Zkj(3) = (xj)/(nj + 2). With these simple considerations it is 
straightforward to obtain results that are close to those obtained with the Sibyll Montecarlo 
code. 

4.2 Comparison with Pythia 

To estimate the systematic uncertainty associated to our calculation we have performed a second 
calculation of the inclusive spectra of unflavored mesons using the Pythia Montecarlo code. Since 
Pythia does not support collisions with a nucleus, we have simulated pp collisions at 10 6 GeV, 
and compared the results with the ones from Sibyll for the same type of interactions. The values 
obtained for the combination Z un fl av (a) defined in Eq. (33) are shown in Fig. |4j It is apparent 



that unflavored meson production in pp interactions is qualititively similar in the Pythia and 
Sibyll codes. Unflavored mesons carry approximately a fraction 0.19 of the projectile energy if 
produced with the Sibyll code; for the Pythia simulation this energy fraction is reduced to 0.16 
(a 20% difference). The energy spectrum in Pythia is however slightly harder. Accordingly, with 
growing a, the Pythia Z-factor decrease a little more slowly, and the difference between the 
models is reduced. For a ~ 2.6 the Z-factors for unflavored meson production calculated with 
the 2 codes coincide. For a ~ 3 the Pythia code gives a result 10% larger. The two codes also 
show some differences in the relative importance of the different mesons. In Pythia the scalar 
(vector) mesons are less (more) important with respect to Sibyll. 

In summary, the description of unflavored meson production in Sibyll and Pythia in pp 
interactions agrees at the level of 10-15%. This level of agreement is however likely to be an 
underestimate of the theoretical uncertainties, because the two Montecarlo codes use very similar 
assumptions. 

The calculation of cosmic ray showers requires also the description of hadronic interactions 
with a nuclear target. The Sibyll code allows to compute the Z factors for both pp and p-air 
interactions. The differences are small but not negligible. For a nuclear target the unflavored 
mesons carry a slightly larger fraction of the energy, but have a softer spectrum (see fig. [4]). The 
first effect is a consequence of the fact that in nuclear interactions the leading baryon is less 
energetic than in pp scattering, and therefore more energy goes into mesons production, on the 
other hand for a nuclear targer all inclusive spectra are softer. 



5 Leptons from charm decay 

Weakly decaying charmed particles (D°, D + , Df, A c and their antiparticles) have a significant 
probability to decay in semileptonic modes such as D° — > K~fj, + v IJj or D° — > K~e + v e , and 
therefore are sources of atmospheric muons and neutrinos. The production of charmed particles, 
however, is dynamically suppressed with respect to the production of pions and kaons, and their 
contribution to the lepton fluxes is usually negligible and remains undetected. On the other 
hand, charmed particles have a lifetime of order r ~ 10~ 12 seconds, and decay with probability 
close to one up to very high energy. The critical energies Ej for D°, D ± , Df, A c are (0.38, 
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0.96, 0.85, 2.4) xlO 8 GeV, several orders of magnitude larger than the ones for pions and kaons. 
Therefore, as the energy grows pion and kaon decay is suppressed and the prompt contribution 
from charm decay will necessarily overtake the standard lepton fluxes. 

The estimate of the lepton fluxes from charm decay has been the subject of many studies 
[El ESI [d HBJ QjJ], with results that span a very broad range. In most cases the prediction 
for this contribution to the muon flux remains always below the one from unflavored meson 
decay discussed in the previous section. If this were the case the charm contribution would only 
be observable in measurements of neutrino fluxes. We do not intend to perform here a new 
calculation of the atmospheric lepton fluxes from charm decay nor a critical review of existing 
predictions. However, we would like to discuss under which conditions this muon flux is above 
the expected one from unflavored mesons. 

Very likely the production of charmed particles does not obey a scaling law of the type 
in Eq. Q. A first order estimate of this contribution can however still be obtained using 
the simple analytic expressions discussed before, treating the hadronic Z-factors as energy- 
dependent quantities. For E% < 10 7 GeV (i.e., a parent charmed particle that decays with 
probability close to 1) and approximating the primary nucleon flux as a power law, the lepton 
flux from charm decay can be estimated as 

"< [t) ~C?(a,E)~ Ma,E)Z ]i (a). (44) 



(K E 



j£{D ,D ,Df,Ac} 



Including the production of charmed particles in nucleon and pion collisions, the factors Aj(a, E) 
are: 

Z Nj (a,E) Z Nn (a) Z vj {a,E) 

A ' ' i-z NN ( a y [i-z NN { a )\[\-z m {a)y y ) 

where for simplicity we have assumed equal cross sections for charm production in p/n or ir^ 
interactions with an air nucleus. 

The contribution from charm produced in pion interactions is suppressed by a factor ZjW; 
and is expected to introduce just a 20-30% correction. Therefore, for a first order estimate 
it is sufficient to model charm production in nucleon interactions. To describe the production 
of the charmed hadron type j in the forward hemisphere we follow the suggestion in [20] and 
parametrize the inclusive spectrum Fnj(x,Eq) as 



F Nj (x, Eo) = ^J>> Pj (nj + !)(!- x)*> , (46) 



This expression can be integrated in the entire interval x € [0, 1]. The corresponding total charm 
cross section is 

ogm ~ A a™(E ) , (47) 

which scales linearly with the mass number A of the target nucleus. The quantity pj is the frac- 
tion of charm events that contain the species j, with YljPj = 1- The -Z-factor that corresponds 



to (46) is 

Z Nj (a, E) ~ -r-j— Pj z{a, nj) - -r-n—z Pj , r~~TTi • ( 48 J 



pA {E) 1J " ' 3 ' o-™\E) " J r(ry + a + l) 



Note that z(l,n) = 1, z(2,n) = l/(n + 2), and z(3,n) = 2/[(n + 2)(n + 3)]. 
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The modeling of the production of leading charmed baryons remains an important open 
problem. This process results into a final state with a A c (the longest lived charmed baryon) 
and, after its decay, high-energy leptons. For this reason we decompose the cross section in 2 
parts, 

a <*- = a A c D + a DD > ( 49 ) 



DD ' 



where the first term accounts for the production of charmed baryons by non perturbative pro- 
cesses. It is natural to expect the energy spectrum of the A c to be significantly harder than the 
spectrum of -D's. Choosing as reference point no ~ 5 and n^ c ~ 1, we obtain 



Cf(a = 3) -1.2 x 1(T 



a vv_ 

DD 



+1.5 x 10 



-6 



100 /ibarn 

a PP_ ■ 

A C D 
100 /ibarn 



z(3,n D ) 
5(3,5) 



+ 



(50) 



z(3,n\ ) 2(3, nn) 

0.6 \,l ^ J +0.4 



5(3,1) 



2(3,5) 



This prediction assumes isospin symmetry, that (primary) scalar and vector charmed mesons 
are produced with the ratio 1/3, and that the production of D s is suppressed by a factor ~ 0.12 
with respect to the production of charmed particles with zero strangeness. 



Comparing Eqs. (50) and (35) one can see that the charm contribution to the atmospheric 
muon flux is below the contribution from unflavored mesons unless the cross section (at Eq ~ 
10 6 GeV) is larger than a^? ~ 100 /xbarn or the energy spectrum of charmed particles is sur- 
prisingly hard. 

The neutrino spectra generated by charmed particle decay can be estimated including the 
appropriate energy spectrum of the neutrinos. Because of the (V — A) properties of the matrix 
element these spectra are a little harder than the corresponding muon spectrum. The decay 
Z-factors (for a = 3) are in the ratios: 



Z Dv ^3) ^ Z DVe (3) ~ 1.25 Z Dtl (3) 
Z KcV (3) ^ Z KcVa {3) ~ 1.16 Z AcM (3) 



(51) 



The different v / 'u ratios for D and A c decay are a consequence of the difference in the available 
phase space in the two cases. In conclusion one expects that, without the inclusion of neutrino 



oscillations, the v e and v„ spectra from charm decay are approximately 



higher that the 



corresponding muon flux. This can be considered a robust prediction, in the sense that it is 
essentially independent of the modelling used for charm production. 

The flux of v T generated by the chain decay of D^r is approximately 30 times smaller than 
the v e or z/„ fluxes. This estimate is obtained assuming cjd s /o"_d — 0.12 and taking into account 
the relevant decay branching ratios and energy spectra. 

The decay probability of charmed particles becomes less than unity for energies larger than 
~ 10 GeV, and therefore the expression in Eq. (44) becomes a poor approximation. To describe 



the resulting lepton fluxes at such energies it is necessary to introduce their critical energies Ej 



and use a functional form of the type given in Eq. (20). 
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6 Photon conversion into fi + fi pairs 



A potentially interesting source of atmospheric muons is the photon conversion into a pair of 
muons: 



(52) 



where Z is the electric charge of an air nucleus. This process is suppressed with respect to the 
production of e + e~ pairs by a factor of {ra e jm^) 2 ~ 2.3 x 10~ 5 . However, high energy showers 
contain a very large number of photons, and it is not immediately obvious that this contribution 
to the muon flux is entirely negligible. 

To estimate the contribution of photon conversions into pairs to the atmospheric 

muon flux it is possible to use the same methods discussed above. The development of an 
electromagnetic shower at high energy is controled by the processes of bremsstrahlung and pair 
production, which are described by the scaling functions 



da, 



e—te-y 



dv 



with v = Ey/E e , and 



da, 



■y—tee 

du 



X — (f(v) 



x A m 



(53) 



(54) 



with u = E e +/Ery. The well known expressions for <p(v) and ip{u) can be found, for example, 
in |22j . In a first approximation the production of muon pairs can be obtained simply rescaling 



Eq. (54): 



da 



7-W 



du 



m. 



m. 



da- 



■y—tee 



du 



(55) 



Assuming a power law for the primary nucleon flux and scaling for the hadronic interactions, 
it is straightforward to obtain that the resulting contribution to the muon flux has the form: 



(K E~ 



(56) 



The a-dependent constant cj?\a) can be calculated as 



C™(a) 



J N^y 



+ 



1 — Znn (1 



Znn)(1 — Znk _ 
A(a — 1) a 7 



A(a - 1) <x 7 - B(a - 1) C(a - 1) 



B(a - 1) 



77 1,. 



(57) 



In this expression one can recognize the effects of an hadronic shower that is the source of an 



electromagnetic one. In Eq. (57) we have left implicit the a dependence of the hadronic Z factors 
Zjy 7 and Z^, the moments of the inclusive photon energy spectra in nucleon and charged pion 
interactions. The photon spectrum is generated by the decay of neutral pions, with smaller 
contributions from the decay of ry mesons and other hadronic resonances. The quantity cj 7 and 
the functions A(s), B(s) and C(s) were introduced by Rossi and Greisen 



(7, 



du ip(u) 



(58) 
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and 



A(s) 
B(s) 
C(s) 



dv <p(v) [1 - (1 
1 

du u s ip(u) 
dv v s ip(v) 



Using in Eq. (57) the numerical hadronic Z-factors obtained with Sibyll at E 
obtain 

cW(a = 3.0) ~ 0.39 x 10~ 6 . 



(59) 
(60) 
(61) 

10 6 GeV we 
(62) 



For a ~ 2.7 the result is approximately 2.5 times larger. 

Comparing this result with Eq. (35) one can see that this contribution to the atmospheric 
muon flux is approximately one order of magnitude smaller than our estimate from unflavored 
meson decay. Most predictions of the muon flux from charm decay are also above the result 
in Eq. (62). Photon conversion into muon pairs is therefore likely to contribute only a small 
fraction of the atmospheric muon flux even at the highest energies. 



7 Summary and Discussion 

In this work we have discussed the main contributions to the atmospheric lepton fluxes at 
very high energy. The standard contribution, due to the decay of charged pions and kaons, is 
suppressed at high energy because most of these long-lived particles interact, and only a small 
fraction (inversely proportional to the energy) decays. This suppression results in very steep 
energy spectra for the lepton fluxes, with a slope approximately one unit larger than the one of 
the primary cosmic ray flux. The zenith angle distributions of these lepton fluxes have also a 
characteristic shape, that for not too large has the well known (cos#) -1 form. Summing over 
particles and anti-particles (and neglecting neutrino oscillations) the lepton fluxes from pion 
and kaon decay are in the ratios: 

/ivi/ M -ri/ e vi/ T ~H 0.35 -f- 0.02 -r (63) 

These ratios reflect the fact that charged pions generate only (fJ>Vu) pairs, with the muon taking 
a larger fraction of the parent particle energy. The flux ratios have an uncertainty significantly 
smaller than the absolute value of the fluxes. The remaining uncertainty is dominated by the 
error in the estimate of the relative importance of kaon and pion production, and is of order of 
15% (20%) for n/vp {vt/vp). 

The contributions to the lepton fluxes from charm decay are characterized by a zenith angle 
distribution that is isotropic for energies below 10 7 GeV, and an energy spectrum that roughly 
follows the primary nucleon spectrum, with corrections related to the energy dependence of the 
charm production cross sections. The prediction has however a large uncertainty due to our 
poor understanding of the dynamics of charmed hadron production. The ratios between the 
lepton fluxes can be predicted with a much smaller uncertainty 

fi -=- -=- v e H- v T ~ 1 -=- 1.2 -=- 1.2 -=- 0.04 (64) 



16 



The decay of unfiavored mesons generates an isotropic muon flux, that follows closely the 
shape of the primary cosmic ray spectrum. Montecarlo codes like Sibyll of Pythia predict that 
this contribution will overtake the standard contribution (for the vertical direction) at an energy 
E ~ 1.5 x 10 6 GeV. Such result reflects that in these Montecarlo codes all mesons that compose 
the scalar and vector SU(3) nonets are produced with probabilites of the same order (taking 
into account a suppression for strange quarks). This implies an abundant production of 77, rf , p° 
and uj mesons, and using the relevant branching ratios one can easily estimate the atmospheric 
muon flux. The prediction can therefore be considered as quite robust. It should be stressed 
that this flux is larger than most predictions for the charm contribution, and therefore it may 
be important for future observations. 

We have also considered the flux generated by photon conversion into muon pairs. The 
muon flux from this source is isotropic and has approximately the same energy spectrum as the 
flux from the decay of unfiavored mesons, but it is approximately 10 times smaller. Photon 
conversion into muon pairs is therefore likely to be of little phenomenological importance. 

It should be stressed that a large uncertainty in the prediction of the lepton fluxes at very 
high energy arises from our poor knowledge of the primary cosmic ray flux. The relevant energy 
range is above the cosmic ray knee at E ~ 3 x 10 6 GeV, where the observations have large errors. 
Still more important is the fact that the composition of the cosmic ray flux in this energy range 
is very poorly known. The prediction of the lepton fluxes depends on the so called nucleon 
spectrum, 

MEo) = J2 A2 MEoA), (65) 

A 

where 4>a(E) is the flux for the nuclear species of mass number A and Eq is the energy per 



nucleon. In Eq. (65 ) the contribution of each nuclear species is weighted by a factor A to account 
for the nucleon multiplicity, and a second factor A is in the jacobian for the transformation from 
total energy to energy per nucleon. If the cosmic ray flux is a power law of exponent a ~ 3 the 
nucleon flux scales with composition oc (A 2 ~ a ) ~ (A^ 1 ). A heavy composition corresponds then 
to a smaller nucleon flux, and to smaller lepton fluxes. 

For the normalization of the nucleon flux shown in Fig. [5] the flux of muons from tt / K and 
unfiavored mesons decay predicted here is of order 



*J and (*U0) ^400 



^^min) * 90 



10 6 GeVl 



Emiri 



10 6 GeV 



(Km 2 yrsr) 1 , (66) 



-2 



(KmVsr)" 1 (67) 



(note the difference in the energy dependence) . For the standard contribution we have performed 
an average over the entire down-going hemisphere. This rough estimates indicate that the 
contribution to the muon flux from unfiavored mesons decay is in principle observable by a 
neutrino telescope of Km 3 . 

The most interesting scientific goal of a large neutrino telescopes such as IceCube is the 
detection of astrophysical neutrinos. The largest signal from astrophysical neutrinos could be 
the isotropic flux from the ensemble of all extragalactic sources. The signatures of such an 
astrophysical neutrino signal are: (a) isotropy; (b) a hard energy spectrum; (c) approximately 
equal fluxes of v e , Vn and v T . The neutrino fluxes generated by charm decay have also the 
properties (a) and (b), and equal fluxes for v e and u^, and therefore constitute a dangerous 
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background. Most theoretical models for the production of astrophysical neutrinos predict an 
energy spectrum harder than what is expected for the charm decay component, however these 
predictions have important uncertainties 

A possible method to separate the astrophysical neutrino signal from the charm decay com- 
ponent is to use measurements of the muon flux to constraint the atmospheric neutrino flux, 
because charm decay generates approximately equal fluxes of muons, and v e . The existence 
of a dominant muon component from unflavored meson decay at very high energy would change 
critically this type of analysis. 
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A Particle Decays 

In this appendix we describe the lepton spectra produced in particle decays that are used in this 
work. The energy spectrum of particle b in the decay of parent particle a can be described in the 
parent rest frame by the function G a b(x), where x = 2E b * /m a and m a is the parent particle mass. 
The spectrum is non vanishing for x in the interval [x m - m ,x mmc \ . The lower limit is x m \ n = 2e 
with e = m b /m a , the upper limit is in general x max < 2 and is determined by the masses of the 
particles in the final state. The normalization of the functions G ab is chosen to be the average 
multiplicity of particle b in the final state, 

(N b )= / dxG ab {x). (68) 

In a frame where a is ultrarelativistic the inclusive spectrum of particle b takes the scaling form: 

^@ a ,E b ) = ±-F ai (y), (69) 

with y = E b /E a and 

F ab (y) = f dx . (70) 

Jy+€ 2 /y Vx 2 -Ae 2 



Equation ( 70 ) assumes that the angular distribution of particle b is isotropic in the parent rest 
frame. This is true for the decay of spin or unpolarized particles. The decay Z-factor is 
defined as: ^ 

Z ab (.a) = [ dy y a - 1 F ab (y) . (71) 
Jo 

In general, particle b can be present in different decay modes. For example r\ decays yield 
muons via the decay modes rj — > n + fi~~f and rj — > /i + fi~ . The functions G ab (x), F ab (y) and 
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Z a b(a) include a sum over all decay channels j where particle b is produced, weighted by the 
appropriate branching ratios Bj. That is: 

F ab (y)= B i F Uv)> Z ab {a)= Y, B i Z lb( a )' (72) 

jgjmodes} j'S{modes} 

with the spectra for each decay mode normalized to unity. 

Two-body decays are very simple to treat. In the rest frame the particle spectrum is 
monochromatic, and in an ultrarelativistic frame it is flat between appropriate kinematical 
limits. For the two-body decay of pions and kaons, such as tt + — > /x+z/^, one has [HE]: 

Fn^y) = 0(y - e 2 ) , F w (y) = [l - % - 1 + e 2 )] , (73) 



with e = m^/m n , and therefore: 



1 - e 2a (1 - e 2 ) 01 - 1 

Z *M = -m 2S > Z -(«) = " " ' ( 74 ) 

a 1 — e z ) a 



For the decays into a pair of particles with the same mass, such as rj — > , one has: 

1 

Vl -4e 2 



F $(V) = r, a o My ~ Vmin(e)] - 9[y - ymax(e)]} , (75) 



with e = m^/m^ and 

ymin,max(e) = 2 (l =F Vl - 4e 2 ) . (76) 

The corresponding Z-factor is: 

Z nu (^) . =R femax 2/min] ■ (77) 

For the decay into three (or more) bodies, the spectrum of the final state particles is deter- 
mined not only by the particle masses, but also by the matrix element of the decay. As a first 
approximation one has that the decay of a parent particle into three massless final state parti- 
cles, just from phase space, has a spectrum (normalized to unity): G(x) = 2x, F(y) = 2(1 — y) 
and Z(a) = 2/(a + a 2 ). To estimate the contribution of unflavored mesons to the muon flux 
one has to consider 3-body decays such as r\ — > /x+/i~7 and uj — > fi + fi~Tt° . The spectrum for the 
latter has been calculated using phase space. For the decay rj — > /i + /i~7, the spectrum without 
matrix element is: 

1 (1 - x) \/x 2 — 4e 2 
C^(ej 1 + e 2 - x 
with Xrr.; n = 2e and x niax = 1 and the normalization factor 



G ^ = 7UT, ^2 - > (^) 



C (e) = ^V'l-4e 2 (l + 2e 2 ) - 6 2 (1 - e 2 ) In 



e 2 (l + v / l-4e 2 ) 



1 - v 7 ! - 4e 2 + e 2 ( v / l -4e 2 - 3) 
In the ultrarelativistic frame the spectrum becomes 



e 2 y 



(l-y)(y-e 2 ) 



(79) 



(80) 
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where the kinematical limits y m \ a and y max are again given by expression ( 76 ) . If the matrix 
element |23| is included then 



\Jx 2 - 4e 2 
(1 - x + e 

with 



G^(x) = A(x)-^——^ + B{x) In 



2e 2 + (1 - x)(x + Vx 2 - 4e 2 ) 
2e 2 + (1 - x)(x - Vx 2 - 4e 2 ) 



(81) 



A(x) = [4- llx + 7x 2 + 2(7-8x)e 2 + 8e 4 ] , (82) 

B(x) = 2 [l - 4e 2 - 2x(l - x)] . (83) 

The normalization C(e) and the distribution in the ultrarelativistic frame FH\P must be evalu- 
ated numerically (see figure [2]). 

To describe the decays of charmed particles we have made the simple but reasonably good 
approximation to treat the dynamics of the charm decay with the matrix element at the the 
quark level for the decay (c — > sl + U() but using the kinematical substitutions: M(c) — > M(D) 
and M(s) -)■ M(K) for D decay, and M(c) M(A C ) and M(s) -)■ M(A) for A c decay. This 
allows to compute the decay spectra analytically. Neglecting the muon mass, in the rest frame 
of the charmed particle the muon and neutrino spectra can be written as: 

, . 1 12x 2 (l-x- e 2 ) 2 , . 

G -wW " W) i-x (84) 

and 

1 2x 2 (l-x-6 2 ) 2 [3-(5-2x)x + (3-x)6 2 ] 

GCU[X) ~ W) {^xf ' (85) 
where x = 2E* V /M(c), e = M(s)/M(c) and 

C(e) = l-8e 2 + 8e 6 -e 8 -12e 4 lne 2 . (86) 



The spectra in a frame where the charmed particle is relativistic can be calculated using ( 70 ) : 
F ceM (y) = { (1 - V ~ e 2 ) [(1 - y) (1 + 2y) - (5 + Ay) e 2 - 2 e 4 ] + 6 e 4 In ( } (87) 
and 

F,M = ^ {5 - - %) - (27 - *V + ^ - ^f^^ 

+6e 4 (3-e 2 )ln 1 --^| . (88) 

An important remark is that the energy spectra of the charged leptons is similar, but slightly 
softer than the corresponding neutrino spectra (fig. [6]). From this one can robustly conclude that 
the expected fluxes of ^ from charm decay will be approximately 15-20% smaller (depending 
on the shape of the spectrum and the composition of the parent charm particles) than the 
corresponding and fluxes. The shape of the neutrino and charged lepton spectra in D and 
A c decay are shown in fig. [6} 

The calculation of the spectrum of the v T produced in the chain decay Dj — > u T + r + — >• 
v T + z7 T + X should include the effects of the polarization of the t + . For the leptonic modes 
(such as t + — > e + v e V T the problem is identical to the case of the chain decay 7r — > jx — > Vn- The 
details for the non leptonic modes will be discussed elsewhere. 
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Figure 1: Zenith angle distribution of very high energy lepton fluxes from pion and kaon decay (thick solid line). 
The dashed line is (cos#) _1 . The thin solid line represents the curve [cos#*(#)] _1 , where 8* (8) is the local zenith 
angle (the angle with the vertical direction) at the point that corresponds to the column density t = 200 (g cm) 2 
for the line of sight defined by zenith angle 8 at sea level. 
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Figure 2: Energy spectra (normalized to unit area) of the muons produced in the 3-body decay of unfiavored 
mesons. The solid lines are for the decay r\ — > /i + /i _ 7 with (thick) and without (thin) matrix element; the dasheed 
line is for the decay u — > n~ir° (using simple phase space). 
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Figure 3: Inclusive spectra of flavorless mesons in p-air interactions (for Eo = 10 15 eV). The 7r° spectrum 
includes the contribution of the decay of unstable resonances. 
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Figure 4: Plot of the quantity Z un n av (a) (defined in equation (33 1) as a function of a. The three curves 
correspond to calculations performed with the Sibyll montecarlo code for pp (dashed line) and p-air interactions 
(solid line), and with the Pythia code (dot-dashed line) for pp interactions. 
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Figure 5: The top line shows the fit of the nucleon flux used in |21) . The curves labeled as fi, and v £ indicate 
the vertical fluxes of atmospheric fx , z/ M + 17^ and i/ e + u e - The calculation of the fluxes extends to higher energy 
of the results of [2T] and coincides with those results for E < 30 TeV. The contribution of the decay of unflavored 
mesons to the muon flux is shown as the thick (red) curve. The points are the measurement of the muon flux of 
the L3 detector fl2l. 
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Figure 6: Energy spectra (normalized to unit area) of the muons and neutrinos produced in charmed hadron 
decays. 
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